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Abstract. A single spherical antenna is capable of measuring the direction and 
' polarization of a gravitational wave. It is possible to solve the inverse problem 

using only linear algebra even in the presence of noise. The simplicity of this 
solution enables one to explore the error on the solution using standard techniques. 
In this paper we derive the error on the direction and polarization measurements 
of a gravitational wave. We show that the solid angle error and the uncertainty 
on the wave amplitude are direction independent. We also discuss the possibility 
of determining the polarization amplitudes with isotropic sensitivity for any given 
gravitational wave source. 
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1. Introduction 



A spherical gravitational wave antenna ideally has equal sensitivity to gravitational 
I waves from all directions and polarizations and is able to determine the directional 

information and tensorial character of a gravitational wave. The solution for the 
inverse problem for a noiseless antenna has been known for some time Q], and an 
analytic solution for an noisy antenna was recently found These solutions are 
quite elegant as they only require linear algebra to estimate the wave direction and 
polarization from the detector outputs. 

By monitoring the five quadrupole modes of an elastic sphere, one has a 
direct measurement of the effective force of a gravitational wave on the sphere 
The standard technique for doing so on resonant detectors is to position resonant 
transducers on the surface of the sphere that strongly couple to the quadrupole modes. 
A number of proposals have been made for the type and positions of the transducers 
^, H]. What all of these proposals have in common is that the outputs of the 
transducers are combined into five "mode channels" gm{t) that are constructed to 
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have a one-to-one correspondence with the quadrupole modes of the sphere and thus 
the spherical ampUtudes of the gravitational wave |0, ^ . 

The mode channels can be collected to form a "detector response" matrix 
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which, in the absence of noise in the detector, is equal to the GW strain tensor, 
expressed in lab frame coordinates. The latter tensor has the canonical form 



(2) 



in the wave frame, and is related to A by an orthogonal transformation — a rotation. 
H clearly has an eigenvector, v^, say, with zero eigenvalue which corresponds to the 
wave propagation direction. The same therefore applies to A, and this enables the 
determination of the wave direction by a straightforward algebraic procedure directly 
from detector data: it is the eigenvector of A with null eigenvalue. 

Things change when the detector is noisy: noise gets added to the signal in the 
mode channels, destroying the equivalence between the data matrix A and the signal 
matrix H. However, it has been shown that under ideal conditions of the noise a 
modified version of the above procedure can be used j^: the eigenvector of the noisy 
A with eigenvalue closest to zero is the best approximation to the actual incidence 
direction of the gravitational wave. 

In this paper we shall be taking an analytic approach to the diagonalization 
problem, whereby errors in the estimated GW parameters can be assessed to any 
desired degree of accuracy. Inherent in this approach is the unambiguous definition 
of the incidence direction estimate, as well as a Cartesian coordinate convention 
for it, which rids us of the ambiguities intrinsically associated to the Euler angle 
characterization for incidence directions near the Poles. Errors in these quantities will 
be shown to be incidence direction independent. In addition, we shall also address 
the problem of estimating the GW amplitudes /i+ and h^, and their errors. Previous 
authors ^ found it impossible to give isotropic estimates of these quantities, a 
very strange result for a spherical detector. We explain why these results come about 
and we show how the problem can be solved by properly including all the necessary 
information. 

The paper is organized as follows. We begin in section |^ by deriving analytic 
expressions for the eigenvalues of A. We then use these expressions to find the 
first order statistical errors in the eigenvalues in section ^, followed by the direction 
estimation error in section ^. Higher order corrections to these errors are presented in 
section In section ^ we discuss the errors on the polarization amplitude estimates. 
We explain why past solutions have direction dependent errors on these quantities 
and we describe a maximum likelihood algorithm, based a hypothesis on the physical 
nature of the source, that fulfills the natural property of source location independence. 
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2. Detector response eigenvalues 

The detector response matrix A is symmetric and traceless and has the eigenvalue 
equation 

Avk = XkVk, k = 1,2,3, (3) 
where the eigenvectors are normalized in the usual way 

Vk-vi = Ski. (4) 
Expanding equation (||) we find 

Xl - g'Xk -D = 0, (5) 
where we have defined 

9^ = 91+92+91+91+91 (6 a) 

D = det(A). (66) 

Solving this cubic equation we find the eigenvalues of A to be 

2 

Afc = --^ 5 cos^fe, fc = 1,2,3, (7) 

where 

+ 2(fc-l)7r 3x/3 D 
Ok = z , and 0086* = -— ^. (8) 

The fc = 3 eigenvalue is identically zero in the absence of noise, so it will generally 
be the one closest to zero in the presence of noise. Random fiuctuations may eventually 
change this (more likely for low SNR), but we shall always take the corresponding 
eigenvector as the best approximation to the direction of the source. The amplitude 
of the wave h can be calculated in many ways from the mode channels (for example, g 
is an estimate for the amplitude), but the semi-difference of the other two eigenvalues 
will give the best estimate Q], 

/i-i(A2-Ai). (9) 

3. Eigenvalue errors 

We assume that the mode channels g„i have uncorrelated noise with zero mean and 
variance (7^^ = E {{Sgm)^}- The lowest order statistical errors in the eigenvalues are 
easily calculated by 

2 ^^fe \^ 2 



where the derivatives of the eigenvalues are given by 
dXk 9m . ^ dD 

T. — = 1= — cos6'fe + - — smOk, (11) 

ogm V3 9 <J9m 

and the derivatives of the determinant D are 
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dii^ 'V^'''' '''''''' ^''"^ 

dD 4 

— = - —g,g, + 2g,g,, (126) 

dD f g5\ , , 

|^ = 2,.(,, + ^)+2,.,3, (12.) 

|£ = (-2g? - 2,92^ + +gl + 2gl) . (12e) 

If the variances on the mode channels are equal then it is easily seen that equations 
(0) and (0) lead to 

crl^^al 1,2,3, (13) 



where is the variance on any one mode channel. Note that from equation (|T3|) all 
three eigenvalues have equal variance to first order. Cross correlations between these 
eigenvalues are also easily calculated, and for equal mode channel variances they are 
equal for all pairs (Afe, A'j.): 

'Jx.y^^-l'Jg, k,k' ^1,2,3. (14) 

Shown in Figs. is this variance as a function of the SNR = S^/^g- Also 
shown is the results of a Monte Carlo type simulation of the errors which take into 
account the higher order perturbations at low values of SNR. As shown, the analytic 
expressions match the simulated errors for high SNR. For low SNR discrepancies arise 
between the analytic and simulated values. This is within expectation since equation 
( p^ is only accurate for large SNR. Higher order correction will be considered below. 

4. Direction estimation error 

We assume that the eigenvector points in the propagation direction of the 
gravitational wave. We want to estimate the fluctuations in the determination of 
this direction caused by the presence of noisy fluctuations in the mode channels. We 
represent with S a difference between a given quantity and its ideal value if there were 
no noise (i.e. Sv^ is the difference between the position calculated from noisy data 
and its real position in the sky). We now take equation (||) for fc = 3 and consider 
fluctuations in it. If these are not too large (high SNR) we can retain only first order 
terms 

[A~\3]6v3^-[5A-S\3]v3. (15) 

This is an equation for Sv^, but the matrix [A — A3] is not invertible. The only 
consequence of this is that we cannot determine the component of Sv^ which is parallel 
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Figure 1. The results of a numerical simulation describing the variance of the 
first eigenvalue for a range of SNR. The dot-dashed line was computed by a f 000 
trial Monte Carlo sirrmlation for the range of SNR. The dashed line is the error 
found from the first order analytic expression and the solid line includes the second 
order corrections. 




Figure 2. The results of a numerical simulation describing the variance of the 
second eigenvalue for a range of SNR. The dot-dashed line was computed by a 
1000 trial Monte Carlo simulation for the range of SNR. The dashed line is the 
error found from the first order analytic expression and the solid line includes the 
second order corrections. 
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Figure 3. The results of a numerical simulation describing the variance of the 
third eigenvalue for a range of SNR. The dot-dashed line was computed by a 
1000 trial Monte Carlo simulation for the range of SNR. The dashed line is the 
error found from the first order analytic expression and the solid line includes the 
second order corrections. 



to ^3 itself. The orthogonal components (those parallel to Vi and V2) can easily be 
found by multiplying equation 15 on the left by Vi and V2 

vi ■ 6V3 — - - — vidAv3, (16) 



V2 ■ 6V3 



Xi — A3 
1 



V1SAV3, 
V2SAV3. 



(17) 



A2 — A3 

An appropriate assessment of the error on a direction measurement is the solid 
angle error AQ. Since |t;3| — 1, this error is given by 

An^nlAvsf. (18) 

where |Ai;3p is the quadratic error in the determination of ^3. To find it we need to 
calculate the expectation of the squared modulus of the above fluctuations, 



|At;3|' ^E{{vi- Svsf + (V2 ■ Sv^f} = E ■ 



V16AV3 



Ai — A3 



V2SAV3 



At — A3 



(19) 



First order calculations only require us to take expectations in SA, while leaving the 
rest untouched, 



6A^^ AmSgr, 



m — l 



where we have defined 



dA 

dgm 



(20) 



(21) 
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Explicitly, 
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We thus have 



5 



9' 



(22&) 



(23) 



Again, setting the variances on the mode channels equal to a^, the sum in equation 
Sj) can easily be done, giving 

2 

|A^3p = 2^. (24) 



From equation (|2J) we see that the error in the incidence direction is independent 
of this direction, as expected of an omnidirectional antenna. Substituting this into 
wc find 

This expression is in perfect agreement with the solid angle estimation error 
found by Zhou and Michelson who used a maximum likelihood technique to estimate 
the wave direction . This is not surprising as the two methods of estimating the 
wave direction have been shown to be equivalent (though the assumptions behind each 
are quite different) 1^. The advantage of our approach is that, by using unit vectors 
(Cartesian components), we are all the time free from the anomalously high errors 
and correlations intrinsically associated to the Euler angle parametrization. 

Shown in Fig. || is the solid angle estimation error as a function of the SNR. Also 
shown is the results of a Monte Carlo type simulation of the errors. As shown, the 
analytic expressions match the simulated errors for high SNR. Deviations however 
appear for lower values of SNR. This again is due to the insufficiency of the first order 
analytical estimates of the errors. In the next section we present improved theoretical 
estimates of the errors by going one order beyond the first in the calculations of 
variances. 



5. Higher order corrections 



In order to improve our theoretical understanding of the error behaviours of the Monte 



Carlo simulations displayed in Figs, 
error terms of equations ( pO| ) and 



we need to go one step beyond the linear 
This requires calculations of higher order 
derivatives for the added terms, and then the resulting general expressions become 
quite complicated. They somewhat simplify for equal mode channel variances, but are 
still rather cumbersome. For example, the next order correction to the eigenvalues. 
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10° 10' 10^ 

Energy Signal-to-Noise Ratio 

Figure 4. The results of a numerical simulation describing the solid angle 
direction estimation error AQ of a source direction measurement due to a finite 
signal-to-noise ratio. The dot-dashed line was computed by a 1000 trial Monte 
Carlo simulation for the range of SNR. The dashed line is the error found from 
the first order analytic expression and the solid line includes the second order 
corrections. 



assuming the mode channel noises are zero-mean independent Gaussian processes, is 
given by (see appendix) 



4 2 

^ '''^ ' ' '''^ '''^ i . = 1,2,3, (26) 



, , dgrn dg„idgm' dg„i' 2 dg,ndg„i' dgmdg^^ ' ' ^ 
After rather long algebra it is found that 



4 



'^l = l<r',-2{l + sin^e,) ^. (27) 

It turns out that the series for cr^^ converges very slowly for low SNR, so that 
equation ( p6| ) only constitutes an improvement on (|l^) for a rather limited range of 
SNR — see Figs. Equation ( p7| ) shows that the errors in Ai and A2 split from the 
error in A3 for low values of SNR, and reproduces the observed behaviour that 
falls below cr^^ and a^^. It is a reasonable approximation to cr^^ for SNR between 30 
and 6, but it is not quite as good as regards cr^ and for those values of SNR. 
Higher order terms would be required for an improvement, but these imply still much 
longer calculation s of derivatives of the eigenvalues up to the fifth order, as can be 
seen in equation ( |A.9|) of the appendix. 



Errors on the inverse problem. 



9 



Similar corrections can be applied to the incidence direction error estimate of 
equation (|2^). They appear to be given by 

|A^3P = 2^ + (28) 

for equal mode channel variances, cr^. Solid angle errors can be directly inferred from 
here: 

^^-JnR i^ + SNR-y\ (29) 

As shown in Fig. ^, the above theoretical prediction is a better approach to the 
behaviour observed in the numerical simulations than is equation (^5|). For SNR less 
than 2, equation ( p9| ) is also insufficient. 

It is important to remind ourselves that for very low SNR the uncertainties on 
the direction estimation are so large that the measurement is almost meaningless. 
Zhou and Michelson decided that a minimum SNR of 10 in energy was necessary for 
a direction measurement Looking at Fig. ^this lower limit seems reasonable, thus 
it is only necessary to have accurate analytical expressions down to that level. 



6. Polarization amplitudes 

We now come to the discussion of errors in the polarization amplitudes h+ and ft-x ■ It 
is easily shown that the uncertainty on the measurement of the polarization amplitudes 
is direction independent if the source position is known ahead of time Q , but it has 
been claimed that in the unknown direction case there is a strong direction dependency 
. This is disturbing given that a spherical detector is equally sensitive to waves of 
all polarization and direction. We argue in this section that the difficulties to find 
isotropic estimates of h+ and hx are ultimately due to the use of unsuitable criteria 
to set up those estimates. We discuss a more natural procedure, based on data set 
processing, that leads to a solution of this problem. 

To understand why a simple estimate of the errors leads to direction dependencies, 
let us look at the basics of the solution to the inverse problem. We are given 3 
eigenvectors and 3 eigenvalues, but these are not independent. The eigenvectors are 
orthogonal, so we actually only get 2 pieces of information from them. We use that 
information to determine the direction of the wave. Next, the strain tensor is traceless 
so the third eigenvalue can be determined from the other two, therefore, we only get 
2 pieces of information from them. We use one of them (actually a combination of 2) 
to get the wave amplitude. Past reasoning suggested that we can use the last piece of 
information to determine the polarization angle. This is wrong. The last eigenvalue 
does not tell us the polarization, but rather is related to the scalar component of the 
GW (or lack thereof). We assume this to be zero for GR, so this can be interpreted 
as a measurement of the non-zeroness of this component. The fact that we decide 
that this should be zero a priori does not give us additional information about the 
polarization, only the level of noise in our system. Without any additional information 
we have an under-determined system which will lead to direction dependent errors as 
seen in reference . 

The additional piece of information needed is the polarization angle a: the angle 
between the GW's axes and the eigenvectors Vi and V2 perpendicular to the incidence 
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direction v^. To this end, we submit our diagonal form of the matrix A to a rotation 
of angle a about to obtain best estimates of /i+ and hy the formulas 

7i_i_ = /i cos 2a, (30a) 
hx = h sin 2a, (306) 

where the pure and uncorrelated noise term in /i+ was dropped out. The GW 
amplitude h = {h\ + h\)^/'^ has been shown to have a best estimate given by equation 
(^. h can be determined with isotropic sensitivity since that is the case with Ai and 
A2, as we have just seen. For a fixed a, equations ( ^OoD and ( |30^ ) give us an estimate 
of hj^ and /ix in the presence of noise in the detector. 
We can use the results of section ^ to see that 

a\ = cos2 2acr^ (31a) 
al = sm^2a(jl (316) 

cr+x = - ^ sin4acr^ (31c) 

and these errors are indeed isotropic, for they only depend on Ug. 

In the absence of further information on the specific physical nature of the source, 
any polarization angle a is valid, for the canonical form of the tensor (||) is invariant 
to rotations about the third axis. A particular choice of a is thus a matter of taste 



in this case, and equations ( pla| )~(31c) give the correct error estimates. A common 
way to resolve the arbitrariness in a is to set the first Euler angle in the rotation 
relating the lab frame to the wave frame equal to zero. 

However, this is a very much observer dependent criterion, for detectors at 
different locations would claim different values for and /ix, even if they agreed 
to be seeing the same source. Errors in /i+ and based on such criterion have been 
shown e.g. by to be strongly direction dependent, which is certainly not surprising. 
It is however paradoxical that a spherical detector should prefer certain directions to 
others to detect a GW signal, therefore this must be reassessed. We now propose a 
more consistent solution. 

It is clear from the above discussion that any criterion to resolve the arbitrariness 
in a, therefore to estimate and ft-x, should be established relatively to the GW 
source, be it known ahead of time or based on a hypothesis to be checked a posteriori. 

Let us, for concreteness, consider a coalescing binary system as the GW source 
[p^ . The signal generated by such a system is given by somewhat complicated 
functions of the space-time variables and a number of system parameters; it will not 
be necessary for our purposes to consider in detail the explicit form of such functions 



(see for example |11 ), it will suffice to use formal expressions indicating the signal 
dependencies: 

h+ = h+{r,t;IC), (32a) 
hx=hx{r,t;IC). (326) 

Here r is the source position, and t is the time. K. stands for the set of 
characteristic source parameters, which in this case include the masses of the stars, 
the inclination of the orbital plane, the semimajor axis, the eccentricity of the orbit. 
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the periastron position, etc. Note that these amphtudes are referred to a set of source 
axes, so they are independent of the detector's location. 

The usual way to estimate the parameters /C is to resort to classical Statistics 
| p2[ , as has been done for example in |l^ for interferometric detectors or in 
for spherical detectors. The fundamental quantity required by such method is the 
likelihood function, A, which is a functional of the (unknown) signal parameters and 
the detector data. 

We then proceed as follows: we construct the likelihood function associated to 



the hypothesis that equations (30q) and (30f) be a fit to equations ( |326| ) and (326) 



for suitable values of the parameters /C. It will thus have the generic form 

A = A(/i;a;/C) (33) 

Standard manipulations of A yield both best estimates of the signal parameters 
/C and of the polarization angle a, as well as errors and cross correlations between any 
pair of these — it is recalled that such are identified as the coefficients of the covariance 
matrix, which is the inverse of the matrix of second derivatives of A [^2|. 

We shall not attempt to give a detailed discussion of this process here. The 
important point to stress is that in the approach just proposed, we have managed 
to have h as the only combination of actual data entering the likelihood function A. 
Errors and cross correlations between parameter estimates will thus ultimately be 
functions only of the errors and cross correlations between the eigenvalue estimates, 
Ai and A2, which we have proved in section ^ to be direction independent. 

So not only a but also the source parameters JC can be determined with isotropic 
sensitivity by means of a spherical GW detector. The same therefore applies to the 
GW amplitudes /i+ and /ix, as indeed expected. 

The quantitative estimation of the errors in h+ and hx cannot however be given 
explicitly until the full parameter estimation problem has been completely solved, as 
interactions between all those estimates will strongly affect one another. 



7. Discussion 



With analytic expressions for the uncertainties on the eigenvectors and the eigenvalues 
of the mode channel matrix (|l|) we can turn our attention to the physical interpretation 
of these values. An unambiguous selection of A3 can be made on the basis that it is 
the third root in equation ([7|). This will usually be the one closest to zero. Then the 
other two represent the amplitude measurement. As proved in [Q, the best estimate 
of the GW amplitude is the semi-difference of these two, (A2 — Ai)/2. 

The third eigenvalue ideally should be zero if general relativity is correct. Once 
noise is introduced this is no longer the case, but the variance on this eigenvalue gives 
us a level of the "non-zeroness" . One can imagine setting a threshold on this eigenvalue 
that is a function of (a function of the SNR) to veto any candidate events that 
have an excessive A3. Many non-GW sources are likely to produce a non-zero A3, 
therefore becoming easily identified and discarded. 

The errors in both eigenvalues and eigenvectors are direction independent. The 
last step in the analysis is the splitting of the GW amplitude h into the usual /i+ and 
hx components. We have shown that this can be accomplished by making suitable 
reference to the source properties, whereby isotropic sensitivity to these quantities 
obtains. This solves the paradox of the anisotropics in the determination of /i+ and 
hx , and stresses the fact that the most fundamental magnitudes to estimate from the 
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detector data are the eigenvalues and eigenvectors of the mode channel matrix: these 
are source independent, and any further progress in signal deconvolution explicitly 
requires reference to the source properties, be them known ahead of time or be them 
stated in the form of a hypothesis to test. 



Acknowledgments 

We are grateful for the hospitality of the ROG group at the Laboratori Nazionali di 
Frascati (Rome), where this work began and partly developed. One of us (SMM) 
thanks E. Mauceli for helpful discussions. Two of us (JAL and MAS) acknowledge 
financial support from the Spanish Ministry of Education, contract number PB96- 
0384, and the Institut d'Estudis Catalans. 

Appendix A. Quadratic error calculations 

Let g = {gi, . . . , gn) be a set of n independent Gaussian random variables, with mean 
/Li — . . . , fin) a-iid variances {i — 1, n). Let f{g) be a regular function of its 
arguments. Because g is a random variable so is f{g), although it will of course be 
generally non-Gaussian. We want to find the mean and variance of f{g), and to this 
end we Taylor expand it around the mean fj,: 

^ fij 5gi Sgj + i 



fig) = fi^J■) + fi % + 77 fij 5gi 6gj + - f,jk Sgi Sgj Sgk + ■ ■ ■ (A.l) 



where 

_ df 



dgidg. 



and 5gi = gi- ^J.i■ (A. 2) 

and the usual convention of summation over repeated indices is adopted in (A.l). 

The mean of f{g) is its expectation value, E{f{g)}, while its variance is the 
difference 

aj=E{f{g)}-[E{fig)}r (A.3) 



Expectation values are to be taken on the basis of the expansion (A.l). Given 
that g is a set of independent Gaussian variables, the expectation of the product of 
an odd number of Sg^s is zero, while 

E{Sg,6g,} = S,, af (A.4) 

E{Sgi Sg-j Sgk Sgi} = 5rj Ski cri o-fe + ^ifc Sji a"^ (t| + Su Sjk (A. 5) 

etc., where no summation over repeated indices is exceptionally assumed in these 
expressions. If the assumption is made that all the g's have equal variances, cr^, then 
one easily finds 

2 fa ^ 8 ^ 48 

where / is a shorthand for /(/i). Likewise, 



E{f{9)} ^f + T^fw^^ + z f^m + 7^ + • ■ • (A.6) 



[E{I{g)}f + ffu <j' + \ {fhjj + 



1 g 

if fiijjkk + 3 fij fijkk) + . . . , (A.7) 



24 
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and 



E{f\9)} = f + iffii + fifi) ^^ + \ iffiijj +^fifijj +^fijfi3) 

1 g 

+ ^ iffiijjkk + 6 fifijjkk + 15 fijfijkk + 10 fijkfijk) O + • • • (A.8) 

So, finally, 



O) = fifi a' + i^^f^J, + - fijfijj a%b 

+ ^ (fifijjkk + 2 fijfijkk + 2 fijkfijkj + ... (A.9) 
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